Residence Time Statistics for Normal and Fractional Diffusion in a Force Field 
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We investigate statistics of occupation times for an over-damped Brownian particle in an external 
force field. A backward Fokker-Planck equation introduced by Majumdar and Comtet describing the 
distribution of occupation times is solved. The solution gives a general relation between occupation 
time statistics and probability currents which are found from solutions of the corresponding problem 
of first passage time. This general relationship between occupation times and first passage times, is 
valid for normal Markovian diffusion and for non-Markovian sub-diffusion, the latter modeled using 
the fractional Fokker-Planck equation. For binding potential fields we find in the long time limit 
ergodic behavior for normal diffusion, while for the fractional framework weak ergodicity breaking 
is found, in agreement with previous results of Bel and Barkai on the continuous time random walk 
on a lattice. For non-binding potential rich physical behaviors are obtained, and classification of 
occupation time statistics is made possible according to whether or not the underlying random walk 
is recurrent and the averaged first return time to the origin is finite. Our work establishes a link 
between fractional calculus and ergodicity breaking. 

PACS numbers: 



I. INTRODUCTION 

Consider the trajectory of a single Brownian parti- 
cle. The total time the particle spends in a given do- 
main is called the residence time or the occupation time. 
The well known example is P. Levy's arcsine law 0, . 
Consider a Brownian motion x(t) = rj(t) where rj(t) is 
Gaussian white noise with zero mean. Levy investigated 
the residence time of the particle in the domain x > 0, 
which we call T + , when the motion is unbounded and 
the total observation time is t. Naive expectation is that 
T + jt = 1/2 with small fluctuations when t — > oo, namely 
the particle occupies the domain x > for half of the time 
of observation. However, instead the probability density 
function of T + jt is given by the well known arcsine law, 
f(T+/t) = [7iy(T+/t)(l-T+/i)]- 1 with < T+/t < 1. 
This probability density has a U shape, which means that 
for a typical realization of the Brownian trajectory, the 
particle spends most of the time in one half of space (say 
x > 0) and not in the other (x < 0). 

Many extensions of this well known result are found 
in the literature. Darling and Kac Im found the lim- 
iting distribution of the time spent in a domain in two 
dimensions, and this line of investigation was extended 
to three dimension by Berezhkovskii et al Q. Lamperti's 
limit theorem gives a very general mathematical foun- 
dation for occupation time statistics (see more details 
in the manuscript). Recently in 7] Pearson's type of 
ballistic motion with random reorientation was consid- 
ered, instead of the usual assumption of an underlying 
continuum process. The basic mathematical theory for 
the calculation of occupation time statistics for Brownian 
motion was developed by Kac, and is usually based on 
the Feynmann-Kac formula (see 0,0 and Ref. therein). 
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Statistics of occupation times is of-course not limited to 
Brownian motion and diffusion, and it is a topic of wide 
investigation , for example in the context of renewal 
processes lEf, theory and experiments of blinking quan- 
tum dots 0, ^| , weak ergodicity breaking of dynamics 
generated using deterministic mapsjl2j , and work distri- 
bution functions of a single spin |l3| . 

The problem of occupation times of a Brownian parti- 
cle in the presence of external field was consider recently, 
by Majumdar and Comtet |l4|. Using the Kac formalism 
PJ they found a backward Fokker-Planck equation whose 
solution yields statistics of occupation times. In the 
problem of occupation time statistics of a particle per- 
forming a random walk on a random walk i.e. the Sinai 
model was investigated. It was shown that statistics of 
occupation times are drastically changed when averages 
over random disorder are made. 

In the first part of this manuscript we consider the 
problem of occupation times for normal Brownian mo- 
tion in an external fields. We solve exactly the backward 
Fokker-Planck equation given in • This solution gives 
a general relation between occupation time and first pas- 
sage time statistics. Besides the theoretical interest in 
such a relation, the solution is used to classify very gen- 
eral behaviors of occupation times based on the corre- 
sponding properties of the first passage times. The later 
are investigated in great detail in the literature Q, and 
we can use this knowledge to solve analytically the prob- 
lem of occupation times at-least for some simple cases. 
For example we show that in the limit of long measure- 
ment times, and for binding force fields, statistics of oc- 
cupation times is determined by Boltzmann's statistics, 
namely the underlying dynamics is ergodic, as expected. 

Statistics of occupation time is important from a fun- 
damental point of view, since if we are able to calcu- 
late statistics of occupation times from some underlying 
dynamics, one can check the validity of the ergodic hy- 
pothesis and its possible extensions. A trivial example 
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is Gaussian Brownian motion in a system of finite size 
< x < L, in the absence of external force fields. Then 
it is easy to show that the residence time in half of the 
system i.e. in the domain (0, L/2), is in statistical sense 
half of the observation time, when the observation time is 
long, as expected. For dynamics described by fractional 
kinetic equations flU ITq . we show that such a simple 
ergodic picture does not hold. 

In the second part of the paper we consider the problem 
of a particle undergoing an anomalous diffusion process. 
We model this behavior using the fractional time Fokker- 
Planck equation [tH Il8| . This fractional framework is 
based on fractional calculus e.g. d 1 ! 2 / 'dt 1 / 2 , which is 
briefly introduced in the manuscript. We show for exam- 
ple, that the general relation between occupation times 
and first passage times we find in the first part of the pa- 
per is still valid, even for the non-Markovian sub-diffusive 
case. Similar to normal diffusion case a classification of 
typical behaviors of occupation times is found, and ana- 
lytical solutions provided. For dynamics in binding force 
fields we find weak ergodicity breaking. In the conclu- 
sions we compare our results on occupation times found 
here using the fractional framework, and recent results 
of Bel and Barkai [T^. I2f1 ] on statistics of residence times 
for continuous time random walks. 

For applications, residence times are of interest in the 
context of chemical reactions 0, Hi| and rather gen- 
erally for statistical analysis of experimental data. Res- 
idence times are very important in the context of single 
molecule dynamics 24, 25]. It is now possible to fol- 
low dynamics of single molecules embedded in condensed 
phase environments, using optical techniques. For exam- 
ple dynamics of single molecules in cells or in solution are 
used to follow chemical reactions in real time, without the 
problem of ensemble averaging found in usual measure- 
ments. A typical experiment uses a laser to investigate 
the dynamics of a particle. In many cases and under cer- 
tain conditions 24] if a particle or a reaction coordinate 
is in a finite domain, the system may emit photons, while 
when the particle is out of the domain the system does 
not emit. Very briefly, the domain width can be imagined 
as the width of the laser beam in single molecule fluores- 
cence experiments when the particle comes in and out of 
resonance with the exciting laser field, due to its diffu- 
sion in space e.g. |2(j, or it could be the Foster radius in 
fluorescence resonance energy transfer measurement e.g. 
[27|. Thus the total time the photons are emitted is ap- 
proximately the residence time, which is proportional to 
the number of emitted photons, which is generally a ran- 
dom variable. For other sources of fluctuations in single 
molecule experiments see |24). 
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where T is the temperature and D is the diffusion coef- 
ficient. As well known the equilibrium of the ensemble 
of particles is the Boltzmann equilibrium, provided that 
the force field is binding. 

Consider a single particle, which at time t = is on 
xq, the observation time of the stochastic dynamics is t. 
The random variable we investigate here is T + , the total 
time the particle occupies the region x > 0. In principle 
during the observation time the particle may cross the 
point x = many times, and then the occupation time 
T + is composed of many sojourn times in x > 0. 

Let P XOi t(T + ) be the probability density function 
(PDF) of T+. The double Laplace transform 



p x ,a (u) 



^-st-uT-' 

e e 



P X0 ,t(T+)dtdT+, (2) 



is defined so that s and t and u and T + are Laplace pairs. 
Majumdar and Comtet [l4| found the equation of motion 
for P Xo ,t(T + ) in double Laplace space 
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F(x ) 



d 2 



Oxq kbT 8xq 



Px ,s (u)—[a + 9 (xq)u] P Xo , a (u) 



(3) 

Where 0(xo) is the step function: 0(xo) = 1 if x > 
otherwise it is zero. This type of equation is called a 
backward Fokker-Planck equation, the operator on the 
left hand side depends on the initial condition xq. Eq. 
© is solved for the matching boundary conditions 



Px ,s(u)\ xo =o+ 



Px o ,s( u )\x o =0- 



dP x „A u ) , 
dx 



=0 + 



8P X0tS (u) , 
dxn 



(4) 



We will re-derive Eq. J3J later as a special limiting case 
of a more general non-Markovian dynamics. 

To prepare for the solution of Eq. JSJ we define the 
following survival probabilities. The probability that a 
particle starting at xq with xq < 0, to remain in the do- 
main x < without leaving it even once, during the time 
t is the survival probability W~ g (t). Let W x ~ (s) be the 
Laplace transform of W~ (t) and similarly the Laplace 
transform of the survival probability in the domain x > 
is W Xo (s) for xq > 0. The key to the solution of the 
problem of occupation times in half space, is to recall 
Tachiya's equation for the survival probability of a par- 
ticle in half space [2II2H 



II. NORMAL DIFFUSION 

Consider a one dimensional over-damped Brownian 
motion in an external force field F(x). The Smoluchowski 
Fokker-Planck equation for the concentration of non in- 



D 



d 2 + FQq) d 
dxZ ksT 8xq 



W~ (s)-sW- 00 = -1 x Q < 0, 

and a similar equation holds for xq > 0. The bound- 
ary conditions for Eq. (J5J are the standard conditions 
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used for the calculation of survival probabilities. Namely, 
120=0 = 0, means that the particle reaches the 



W- [s) 



boundary on x = instantaneously if the particle starts 
very close to the absorbing boundary and if xq — > — oo 
survival is unity. 

The solution of Eq. © for the occupation time is 

P xo , s («) = W-(s) + [1 - sW-(s)] G s (u) 

if xq < 0, 

P XQ , S (u) = W+(s + u) + [1 - (s + u)W+(s + u)] G s (u) 

(6) 

if xo > 0. From Eq. (jfjj) the physical meaning of 
G s (u) becomes clear, it is the double Laplace trans- 
form of Gt(T + ) the PDF of the random variable T + 
for a particle starting on x = 0. The PDF G t (T+) 
contains the information on the problem of occupation 
times, while the survival probability was investigated pre- 
viously by many authors, hence in what follows we in- 
vestigate Gt(T + ). Using Eq. the reader can verify 
that Eq. © is indeed the general solution of the prob- 
lem of occupation times Eq. Using the boundary 
condition W~ Q (s) = W x (s) — and the solution 
Eq. JBJ it is easy to see that the boundary condition 

p x ,s( u )\x =0+ = P xoA u )\x o =0- = G s (u) in Eq. gj is 
satisfied. 

The second matching boundary condition in Eq. igj, 
on the derivatives of P XQyS (u) yields G s (u) using Eq. © 



G s (u) = 



J+(s + u) - J~(s) 
(s + u)J+(s + u)- sJ-(s) 



(J) 



where the currents are 
J+( s+u ) = J -\ xo=0+ , J-(s)~ 



dx 



dx 



(8) 

Eqs. 10 CJ) are the main results so far since they yield 
a general relation between statistics of occupation times 
and survival probability currents. From Eq. ||HJ) we see 
that the solution of the problem of occupation times is 
found in terms of two solutions of the corresponding first 
passage time problems, the first for a particle starting on 
x > and absorbed oni = (i.e. J+) and the second 
for a particle starting on xo < and absorbed on x = 
(i.e. J~). Thus the problem of residence times is solved 
in three steps: 

(i) Find solutions of two first passage time problems for, 
xq > and xo < in Laplace space. 

(ii) Use Eq. Q) to find the solution of the problem of 
residence times in double Laplace space. 

(iii) And then use a two dimensional inverse Laplace 
transform to get Gt(T + ) from G s (u). 

Since there exists a vast literature on the solutions of the 
problem of first passage time , the relationship Eq. Q) 
is very useful for the calculation of statistics of occupa- 
tion times. We note that some connections between first 
passage times and occupation times, which are different 
and in our opinion less general than Eq. JSJl, appeared 



previously in the literature 5, 8, 20]. Finally, while we 
considered the occupation time in half space, occupation 
times in a finite domain are also obtained in a similar 
way, and it is straight forward to extended our results to 
higher dimensions. 

Majumdar and Comtet classify statistics of occu- 
pation times according to behavior of the potential field, 
in particular they consider motion in stable, unstable and 
flat potential fields. Here the relation between survival 
currents and statistics of occupations times, Eq. Q) can 
be used to characterized certain very general and new 
behaviors of occupation times. 

Survival probabilities in a finite and infinite domain 
exhibit three well known typical physical behaviors 0], 
we consider the right random walk (i.e. xq > 0) and simi- 
lar classification holds for the left random walk. Later we 
will classify behaviors of residence times based on these 
three behaviors of first passage times. 
Case 1 The random walk is recurrent, and the average 
first passage time from Xo to is finite. Such cases cor- 
respond to diffusion in a system of finite size, when the 
particle cannot escape to infinity, e.g. the driving force 
field is binding. 

Case 2 The random walk is transient, i.e. the survival 
probability in x > is finite in the limit of long times. 
Such cases happen when the external force drives the 
particle far from the origin, and the system is infinite. In 
that case in the limit of small s 



w: 



(9) 



where ej is the survival probability of the particle start- 
ing on xo, without reaching x = 0, when t — > oo. Similar 
notation is used for the left random walk, with e~ . 
Case 3 Random walks are recurrent, though the average 
first passage time is infinite. A particularly common situ- 
ation is the case when the survival probability decays like 
t -1 / 2 for long times. This happens if the non-diverging 
external field F(x) — for x > x c and the system is in- 
finite, namely when diffusion controls the long time dy- 
namics. For such a case 



w: 



A t (> 

„l/2 < 



0. 



(10) 



where A+ Q > depends of-course on the details of the 
force field. 

We now consider certain general properties of the statis- 
tics of occupation times for the three cases. 

Case 1 The long time behavior of Gt(T + ) is now con- 
sidered. When the left and right random walks, starting 
at xo < or xq > 0, respectively, are recurrent and the 
average first passage time is finite. For this case the small 
s limit yields 



W x k o (s = 0) = (4), 



(11) 



where (t x ) is the average time for the particle starting 
on xq < (or xq > 0) to reach the origin for the first 
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time. The small s and u limit, their ratio arbitrary, of 
Eq. Q gives the long t and T+ behavior of G t (T + ), we 
find 



G s (u) 



1 



8(tZ > 
x ^0 ' I 
Sin 1^0 = ° 



(12) 



9 <*5o> l „ 

Sin 1^0 = ° 



8 <*»0> | „ 



The differential equation for (ti.) is well known and is 
obtained from the small s expansion of Eq. JSJ) 



D 



dx 2 



£(xo)_a_ 



= -1. 



(13) 



Solving this equation, using a similar equation for (t~ ), 
and inverting Eq. Ijl2|l to the time domain we find the 
expected ergodic behavior 



5 (T 



(14) 



where Pit is Boltzmann's probability of occupying x > 



Jo 



' da; 



Z 



(15) 



Z = exp — da; is the normalizing partition func- 
tion and U (x) is the binding potential, with F(x) — 
-dU(x)/dx. 

Case 2 We consider the case where both the left and 
the right random walks are non recurrent. The survival 
probabilities in the two domains are e+ and e~ , in the 
long time limit. Then using Eqs. I|7I9|I we find for t — > oo 



G t (T+) 



a 



-5{T + )+a + 6(T+ -t) 



where 



9x |a: ° =0 



de ig I _ de *o I 

9x l x o=0 axo |x =0 



(16) 



(17) 



and a - = 1 — a + . Since the particle always manages 
to escape either to the left or to the right, eventually 
the particle will either reside in the left domain or the 
right domain forever, hence the delta functions in Eq. 
(jTBj) . The weights of these delta functions are given by 
the derivatives of the survival probabilities only. 

Case 3 We now consider a case where both the left 
and the right random walks are recurrent, though the 
average first return time from xq to x — is infinite, in 
such a way that Eq. Q10(l is valid. Then in the small s 
and u limit 



G s {u) 



-1/2 



• n(s- 



-1/2 



S V» + ft(« + ti)V2 



where the asymmetry parameter is 

9A+ , 



TZ= - 



dxp \xo=0 
dx 1^0=0 



(18) 



(19) 



Transforming to the time domain we find the asymmetric 
arcsine PDF 



where 



/(*) 



1 /T+ 



K 



* x 1 / 2 (1 - x) 1 ' 2 [TZ 2 (l-x)+: 



(20) 



(21) 



When 1Z — 1 we find the arcsine law. Note that the PDF 
Eq. diverges on T+/t = 1 and T+/t = 0, hence 

events where the particle always occupies (or hardly 
never occupies) the domain x > have a significant con- 
tribution. 

Another general result obtained from the small u ex- 
pansion of Eq. J2J is for the average occupation time 



J+(s) 



s 2 [J+( s )-J-(s)} 



(22) 



where C~^ t is the inverse Laplace transform. If the po- 
tential field is binding and the random walk is recurrent, 
then (T + ) ~ P^t- Similar relations between higher or- 
der moments of the occupation times and the survival 
probabilities are obtained in a similar way. In Sec. IIVI 
we consider several particular examples, which explain 
in greater detail the meaning of the general results ob- 
tained in this section. First we generalize our results for 
fractional dynamics. 



III. ANOMALOUS DIFFUSION 

Anomalous diffusion and relaxation is modeled based 
on the fractional time Fokker-Planck equation (FFPE) 
[ItI Il8[ . the concentration of non-interacting particles 
obeys 



d a c(x,t) 
dt a 



= D a 



dx 2 dx k b T 



d F(x) 



c(x,t), 



(23) 



where D a is a generalized diffusion coefficient and < 
a < 1. A brief mathematical introduction to the FFPE 
is given in Appendix A. We recall physical properties 
of the FFPE. (i) when F(x) — and for free bound- 
ary conditions we have the fractional diffusion equation 
[30T l3lL 13^ . I33L l34l | with anomalous diffusion (a; 2 ) cx 
t a . (ii) In the presence of a binding time independent 
force field the equilibrium is the Boltzmann distribution 
[rH Hsj . (iii) Generalized Einstein relations are satis- 
fied in consistency with linear response theory 0, . 

(iv) Relaxation of modes follows the Mittag Leffler de- 
cay, related for example to Cole-Cole relaxation 0, 0] . 

(v) In the limit a — * 1 we recover the standard Smolu- 
chowski Fokker-Planck equation. The FFPE is derived 
from the continuous time random walk [Tsj . Its mathe- 
matical foundation is P. Levy's generalized central limit 
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theorem applied to the number of steps in the underly- 
ing random walk [3^ |3(j . A very general solution of the 
FFPE in terms of the solution of the standard a = 1 
Fokker-Planck equation was given in [3S| (i.e., subordi- 
nation, and the inverse Levy transform) . Recently there 
is some controversy on how to apply boundary conditions 
[37L IssL l8?il] for the anomalous case. Applications of frac- 
tional diffusion modeling include: Scher-Montroll time of 
flight transport of charge carriers in disordered medium 
|3fij . dynamics of ion channels (^3, relaxation processes 
in proteins j^JJ , and dielectric relaxation 1 an d deter- 
ministic diffusion |32| . For a review and a popular article 
on fractional kinetics see [TH Hr| . 

Similar to the normal diffusion case, we define 
P X0 .t(T+) as the PDF of T+ and P xo , s (u), its double 
Laplace transform. As we show in next subsection the 
differential Eq. for P a 
by the FFPE Eq. @ 



x ,s(u) for the dynamics described 
is 



D a 



d 2 F{x Q ) d 



Ox 2 



p x ,s (u) 



[s + e (a*) u] a p Xo , s (u) = - [s + e (ar ) ■ (24) 

The boundary conditions for Eq. (|24[1 are identical to the 
normal diffusion case a — 1 given in Eq. Q}. Eq. I|24l) is 
important since as pointed out in [19| , fractional dynam- 
ics is weakly non-ergodic [4^ |. namely occupation time 
statistics is not described by Boltzmann equilibrium even 
in the limit of long time and for binding potential fields. 
Thus the FFPE l|23|) cannot be used to describe time av- 
erages of physical observable due to ergodicity breaking, 
and the interpretation of results derived from the FFPE 
must be treated with care. Eq. I|24ll is a remedy for this 
problem since as we will show it yields a fractional frame- 
work for the calculation of non-trivial distribution of oc- 
cupation times [i.e. generalization of Boltzmann's statis- 
tics Eq. dJ]. Eq. is a fractional backward Fokker- 
Planck equation in double Laplace space, formally one 
may invert it to the time domain using material frac- 
tional derivatives 0], however in practice we solve this 
equation in double Laplace space and only then make the 
inverse double Laplace transform. 

Interestingly the solution of the fractional Eq. I|24[l is 
identical to that found for normal diffusion case, namely 
our main results Eqs. (|6I7I8|I arc valid also in the non- 
Mar kovian domain < a < 1. Now W^r (s ) needed for 
the calculation of J (s), is the Laplace transform of the 
survival probability for the fractional particle. Thus Eqs. 
(|6I7I8|I have some general validity beyond normal Marko- 
vian diffusion. 

To prove that Eqs. [KilSfl are still valid we must first 
find the differential equation for W^(s): the survival 
probability of a fractional particle starting on xq > 
in the domain x > 0. We can prove that 



D a 



dxl 



F(x ) d 
ksT dxo 



W+ (s) - s a W+ (s) = - 



and a similar equation holds for xq < 0. Eq. I|25|) is 
the fractional generalization of Tachiya's Eq. (JSJ. The 
derivation of Eq. (|25|l is based on results obtained in 
|35| and is simple once the sub-ordination trick is used 
(see some details in Appendix A). Now using Eq. 125|) 
it is easy to verify that Eqs. KilSI) are solutions of the 
fractional Eq. (|231 . 



A. Derivation of Fractional Equation for 
Occupation Times 

In this subsection we derive our main result Eq. I|24|) 
using the assumption that the underlying dynamics is 
described by the fractional Fokker-Planck equation ill' Ml) . 
The latter describes long time behavior of the continu- 
ous time random walk (CTRW), which is the underlying 
random walk process we have in mind. In the CTRW a 
particle performs a one dimensional random walk on a 
lattice, with jumps to nearest neighbors only and with 
random waiting times between jumps. In CTRW the 
waiting times between jumps are independent identically 
distributed random variables, namely the CTRW process 
is renewed after each jump. The PDF of waiting times 
is ip(t). Two classes of CTRWs are usually considered, 
the case when the average waiting time is finite, and the 
case when ip(t) oc when t — > 00 and < a < 1. 

The latter case leads to a non-stationary behavior, a ging , 
anomalous diffusion and weak ergodicity breaking [43j. 
The lattice spacing is e, eventually we will consider the 
continuum limit when e is small. On each lattice point 
we assign a probability for jumping left and a probabil- 
ity of jumping right. This dynamics in the continuum 
limit leads to behavior described by the FFPE, when de- 
tailed balance conditions are applied on the probabilities 
for jumping left or right |l8j (i.e. probabilities to jump 
left or right are related to external force field and tem- 
perature). In what follows we start with some general 
arguments, assuming only a renewal property of the ran- 
dom walk, without limiting our selves to a specific model. 

The random position of the particle is x{t). The total 
time the particle spends on x > is T + , i.e. the occu- 
pation time of half space. The particle starts on xq and 
assume that x$ > 0, later we generalize our results to the 
case xq < 0. We define the PDF of first passage times, 
from xq to x — — e, as tp^ (t). The PDF of first passage 
times from x = to x = — e (x = — e to x = 0) is denoted 
with ip + (t) and \^~{t)\ respectively. 

We assume that the first passage times PDFs ip + (t) 
and tp~ (t) do not depend on £0 and that sojourn times 
in domain x > and x < are statistically indepen- 
dent. Such assumption holds for Markovian dynamics 
but is not obvious otherwise. For CTRW dynamics the 
assumption is correct, since as mentioned the CTRW pro- 
cesses is a renewal process. The process is mapped on a 
two state process 



(25) 



%{t) = 



1 x(t) > 
x(t) < 



(26) 
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and hence T+ = f Q 6 x (t)dt. Since either the particle is 
in the domain x < or not the dynamics is described by 
a set of sojourn times 

Here the PDF of t Xq is V^oW, the PDF 01 is V>~(*)i 
the PDF of is ip + (t), etc- All the sojourn times are 
assumed mutually independent, which means the process 



is renewed once the particle jumps from x = 0tox = — e 
or vice versa. 

Let f Xo , t (T+) be the PDF of T+ when the total ob- 
servation time is t. Let f Xo ,s(u) be the double Laplace 
transform of f XOy t(T + ). A calculation, using methods of 
renewal theory, similar to the work of Godreche and Luck 
yields 



fx ,s(u) 



l-iP+(s + u) 



s + u 



■ l-i/j+js + u) l-jy-js) 

V (s) : 1 



s + u 



1 



1 - ip+(s + u)ijj-(s) ' 



(27) 



Where -0+ (s + u) = J Q exp[— (s + u)t]^+ (t)dt is the Laplace transform. If i/>+ (t) = ?A + (t) = ip (t) we recover a result 
in 9J. If the particle starts on xq < then one can show 



fx A u ) = + 



t/> + (s + u) 



l-if}~(s) l-ip + (s + u) 



1 - ip+(s + u)ip-(s) 



(28) 



We now consider the case when the underlying dynamics is described by the FFPE. By definition the first passage 
time PDFs are related to survival probabilities according to 



W+(t) = l- / 7/>+(f)dt, 

Jo 



or using the convolution theorem in Laplace space 

Hence we may rewrite Eq. Q27JI 

/*„,.(«) = W+(s + «)+[!- (s + u)W+(s + u)} U-(s) 



l-^+(s + tt) 1-V>~(s) 



1 



1 - '0+(s + u)lp-(s) 



(29) 



(30) 



(31) 



Notice that f XOtS (u) depends on xo only through the sur- 
vival probability W XQ (s + u). If we apply the backward 
Fokker-Planck operator 



D a 



d 2 | F(x ) d 
dxZ ksT 8xq 



on this equation and use Eq. J2HJ, namely we assume 
that the underlying dynamics is described by the FFPE 
in the continuum limit, we obtain at once our main result 
Eq. H24f> . Similar method is used for the case xo < to 
complete the proof. 

We now derive our main Eq. Q using the continuum 
approximation. We consider the case when the particle 
start on Xq — hence we have ipj = ?/> + and define 



G s (u) = f Xo =o, s (u). 



(32) 



Generally G s (u) is not identical to G s (u) and our aim 
now is to find the conditions when these two functions 



are identical. Using (|27|) we have 

'l - i'+(s + u) 



G s {u) 



l-^-(s) 



(s + u) 



l—ip+(s + u)tJj-(s) ' 



(33) 



In the continuum limit we have the following e expansion 



Ipe (s) - ^ e= o( s ) 



<hf£ 

dx 



(34) 



and a similar expansion holds for ^+(s). Where the sub- 
script e in Eq. is added to emphasize that the PDF of 
the first passage time from lattice point — e to the origin 
0. Note that ip7=o( s ) = 1> smce the particle on the origin 
is immediately absorbed. Inserting the expansion l|34|) in 
Eq. (|33|l and using a similar expansion for + u) we 
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find that when e — > 

1 dji + (s+u) 



G 8 (u) 



1 di>-(s) I 

Js+u) dx \xo=Q s dx Q 1^0=0 

9l/l+(s+«) 1 (s)| 

|x =0 gj^ |x =0 



(35) 



9s 



Using Eqs. I|3UI35[) we obtain our main result Eq. (JJJ 
and G s (it) is identical to G s (u) in the continuum limit of 
e -► 0. 



IV. EXAMPLES 
A. Weak Ergodicity Breaking 

We now consider anomalous dynamics in a binding po- 
tential field U(x), e.g. U(x) = kx 2 /2 with k > 0, and 
F(z) = — dU(x)/dx. As we showed already the long 
time behavior of Gt(T + ), for normal diffusion, yields an 
ergodic behavior Eq. (|14fl . Residence time statistics for 
sub-diffusion in a binding potential field was considered 
previously in [l9l | , using the continuous time random walk 
model on a lattice. Here we consider a fractional Fokker- 
Planck equation approach showing that concepts of weak 
ergodicity breaking in [19J are valid also within the frac- 
tional framework. 

In |35| it was shown that the random walks in a binding 
field are recurrent for a < 1, just like the normal case 
a = 1. Using the subordination trick (see Appendix), or 
analyzing Eq. (|25|l we find that for s — ► 



(36) 



When « = lwe have (t^)" = (t^f ) and the behavior in 
Eq. pip. For < a < 1 we have in the time domain 
W^ (t) oc t~ a , reflecting the long tailed trapping times of 
the underlying CTRW. The (r^r ) " are amplitudes which 
satisfy 



D a 



o 2 



F(x ) d 



dx 



ksT 8xq 



-i. 



(37) 



This equation is obtained from the small s expansion of 
Eq. l|25Jl. Using Eq. and Eq. (J7J we find in the limit 
of small s and u 



G s {u) 



K(s + u) c 



K(s + u) a + s a 
where the asymmetry parameter is 



11 = 



dx Q 1 2:0=0+ 



aor- 



(38) 



(39) 



dx \x =0 



Inverting to the time domain, we see that the PDF of T + 
in the long time t limit is described by Lamperti's limit 
theorem lo| 



(40) 



where the scaling function is 

S a (TZ,p) = 



sm ira 



Up 01 - 1 (1 -p) c 



Ti" K 2 (l-p) 2a +p 2a + 211(1- p) a p a cos na 

(41) 

This function is normalized according to 
f Q 5 a (K,p)dp = 1. When a = 1 we find the er- 
godic behavior in Eq. (|14|) , while clearly if a < 1 we find 
a non-ergodic behavior. The parameter 1Z is called the 
asymmetry parameter. It can be calculated solving Eq. 
$1). we find 



^-L=o + = 4- f e -[W)]/^, (42) 



- L /" p -[tf(i')-£f(0)]/Wjj 



(I 



Using Eqs. (|39l42l43|) we find 

Pi 



K 



1 -pi 



(43) 



(44) 



where Pg is Boltzmann's probability of finding the par- 
ticle in the domain x > Eq. {T5J. Eqs. (|40I44|I were 
found previously in [Tflj using a different approach. One 
can show that the average occupation time is 

(T+) ~ P+t (45) 

and fluctuation are very large if a < 1 

(T+ 2 > - (T+> 2 ~ (1 - a) P+ (1 - P+) t 2 . (46) 

We inrtoduce a measure for ergodicity breaking the EB 
parameter 



(T+ 2 ) - (T+) 2 1 
VB= { ^ ' ~(l-a)- 



- P+ 

B 

P+ ' 

B 



(47) 



which is zero in the ergodic phase a = 1. 



B. Diffusion in an interval 

We consider the case where the particle is free to dif- 
fuse in an interval of total length L + + L~ . The particle 
is initially on the origin x — and reflecting boundary 
conditions are on x = L + and x = —L~. Statistical 
properties of T + the time spent in (0, L + ) are now inves- 
tigated. 

The survival probability can be calculated using Eq. 

M 



w + (t 

" xo y 



^ cosh [y^$|^(Z, + — so)] 
cosh ^ 75— L + ) 



(48) 
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and when a — 1 we recover a text book result 
Eq. we find 

G s (u) = 



Using 



(s + u) a ^ 2 1 tanh 


'(s+u) a/2 L+~ 


+ s Q / 2 - 1 tanh( s ^i 


--) 


(s + u) a ^ 2 tanh 


' (s+u)"/ 2 L+~ 


+ W 2 tanh 







(49) 

For free boundary conditions, namely in the limit 
where the system size is infinite L + — > oo and L~ — > oo 
we find 



G fl (u) 



(s + u) Q/2 + W 2 



(50) 



Inverting to the time domain, the PDF of T + is the sym- 
metric Lamperti PDF with index a/2 



G t (T+) = - t 5 a/2 (l, 



(51) 



When a = 1, i.e. the case of normal Gaussian diffusion, 
we recover the well known arcsine distribution. As a is 
decreased we are more likely to find the particle localized 
in x > or x < for a time of the order of the obser- 
vation time. Indeed when a -* the PDF of T + is a 
combination of two delta functions with T + = t (particle 
always on x > 0) or T + = (particle always in x < 0). 

A different behavior is found for finite L + and L~ , then 
using Eq. (|4*9")l we find 



1 / T+ 
G t (T+)~- t 6 a/2 (l,— 



t « 



mm(L+,L-) 2 



A, 



l/a 



(52) 



For these time scales the particle does not interact with 
the boundaries, and Gt{T + ) is the symmetric Lamperti 
PDF with index a/2. In the long time limit, correspond- 
ing to small s u limit we find using Eq. H49|) 



G S {U) r 

and hence when t 



(s + u) a - l L^ 



(s + u) a L+ + s a L- 



oo 



(53) 



(54) 



This is in agreement with our more general result Eqs. 
(|40I44|I namely for the case of free diffusion Pg = 
L+/(L+ + L ) and hence TZ = L+/L~. If L+ ^ L~ 
the PDF of T + is as expected non-symmetric, reflecting 
the tendency of the particle to reside in the larger inter- 
val [say (0, L + ) if L + > L~] for longer times compared 
with the shorter domain. For long times an equilibrium 
is obtained: for a — 1 an ergodic phase is found where 
T + jt = L + /(L~ + L + ) while for a < 1 weak ergodicity 



breaking is found. We see that the statistics of occupa- 
tion times exhibits a transition from a symmetric Lam- 
perti PDF with index a/2 when diffusion is dominating 
the dynamics, i.e. for short times, to a generally non- 
symmetric Lamperti PDF with index a, for long times 
when the particle interacts with the boundaries. Such a 
transition is not limited to free diffusion as we show later. 

If L~ — > oo while L + remains finite a different behav- 
ior is found. Now the particle can be found either in a 
domain of finite length < x < L + or in the infinite 
domain — oo < x < 0. Statistically we expect of-course 
that the particle will reside more in i < 0, though the 
random walk is recurrent hence, after each sojourn time 
in x < the particle is ejected back to x > 0, provided 
that we wait long enough. However the average return 
time from a point in x < to some point in x > is 
infinite, and this means that simple scaling T + ~ t does 
not not hold. For this case we have 



G s {u) 



(s + u) a ^ 2 1 tanh 


"(s+u) Q/2 L + " 


+ W 2 " 1 






(s + u) a ^ 2 tanh 




+ W 2 



To investigate deviations from simple scaling we consider 
moments of the random variable T + , using the small u 
expansion of Eq. I|55|l . The average occupation time in 
< x < L + is 



(T + ) = C- 



1 

2? 



1 - e 



(56) 



where C~l+ t is the inverse Laplace transform. This ex- 
pression is inverted using one sided Levy stable functions, 
recall 



a/2,2L+ /VA 



■,l(*)=^ 



-2 S °/ 2 l+/ v / d; 



(57) 



and see |35| and Ref. therein for more mathematical 
details on this function. The one sided stable cumulative 
distribution is 



' -JdZ Jo ^/ol 



and hence 



L a /2, 2SL,l (*) 



dt. 



(58) 



(59) 



For short times (T + ) = t/2, since then the particle does 
not have time to interact with the boundary, and it 
spends half of the time in x > 0. For long times 



(T 4 



sfDZT(2-a/2Y 



(60) 



We see that as the process becomes slower, namely when 
a is decreased, the particle tends to stay more in < 
x < L+ i.e. (T+) oc t for a = but (T+) oc t 1 ' 2 if 
a = 1. We explain this result for normal diffusion a = 1 
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by thinking about the process as a two state process, i.e. 
the particle is either in x < or in x > 0. Sojourn times 
in x > are finite, since the interval < x < L + is 
finite. The PDF of times in state x < follow the t -3 / 2 
power law tail due to usual diffusion. The number of 
time the particle will cross zero is (n(t)) ~ t 1 / 2 at-least 
for a lattice CTRW process (in the continuum limit this 
question is not well defined). Hence we expect (T + ) = 
(n(t)) * average time in < x < L + ~ t 1 ! 2 as we find. If 
a < 1 we expect (n(t)) ~ t a ^ 2 and the average time is < 
x < L+ is proportional to f tp(t)tdt ~ /' H 1+a )*di ~ 
t 1 ^", hence we get (T + ) oc i 1- "/ 2 . The second moment 
of T+ is 



(T+ 2 ) -4(1 -a) 



L+ t 2 - a / 2 
/D^T(3- a/2)' 



(61) 



Note that we do not have simple scaling and (T + 2 ) cx 
t 2 ~ a / 2 is not proportional to (T + ) 2 cx t 2 ~ a , and hence 
the PDF of T + does not have a simple scaling. 



namely a symmetric Lamperti PDF with index a/2 de- 
scribes the residence times. Such behavior is independent 
of the drift and can be understood if we notice that for 
short times the dynamics is governed by diffusion not 
drift. To see this recall that the scaling of these two pro- 
cesses is x — t a / 2 (diffusion) and x — t a (drift) and hence 
for short times diffusion wins. For long times we use the 
small s , u expansion of Eq. (|65|) , G s (u) — l/(s + u) which 
gives the expected behavior G t (T+) — 5{T + — t). 
The mean occupation time is 



{T + )=C-i t 



11 + VI + is a T a 



, 2s 2 VI + 4s q t q 
For the case of normal diffusion a = 1 we find 



The long time behavior is 



- - 1 ) W- 

2r J V 4r 



(T+> ~ t 



(68) 



(69) 



(70) 



C. Diffusion with Drift 

We consider anomalous diffusion in the presence of a 
constant driving force F > 0, for an infinite system. The 
biased diffusion yields a net drift (x) = D a Ft a /[k b TT(l+ 
a)]. Since F > the particle will escape to infinity, hence 
for a particle starting on x — we expect T + — t when 
t is large. 

The survival probability in the right half space 



W + (i 



1 — exp 



Ff3+(s)x 
2k b T 



(62) 



when xq > 0. To obtain the survival probability for left 
random walks replace (3 + (s) in Eq. I|62|l with /3 _ (s), and 



where 



Using Eq. (JSJ) 



P ± (s) = 1 ± Vl + 4s Q r c 



(k b Tf 



F 2 D a 



G s (u) 



(s+u.) s 



[3+(s + u)- [3-(s) 
For s a r a » 1 and u a r a » 1 

(s + u)"/ 2 - 1 + W 2 " 1 



(s + u) a > * 
Thus for short times t «t 



a/2 



G t[ T+)^U a/2 (l T -) 



(63) 



(64) 



(65) 



(66) 



(67) 



the leading term (T + ) — t is expected since as mentioned 
for long times the particle is always in x > when F > 0. 
For short times 



1 



1/2 



1 



3/2 



■0(f 



5/2N 



3^7? \T/ 30-^7? \T 

(71) 

The leading term (T + ) — | shows that at short time 
diffusion not drift is dominating the process, hence from 
symmetry half of the time the particle is on 1 > 0. For 
the sub-diffusive case a < 1 we investigate the long time 
behavior of (T + ) using the small s expansion of (|68|) and 
then inverting to the time domain 



t 



1 - 



1 



(k b T) 2 



(72) 



r (2 - a) F 2 D a t a 

At short times we use Eq. (|68|) and Hankel's contour 
integral, for the T(z) function, and find 



1 



1 



FjD r 



2T (2 + a/2) ( feh T) ; 



a/2 



(73) 



We note that results for the case F < can be easily 
obtained from our results for F > 0. The distribution 
of times T~ in x < when i 7, > is equal of-course to 
the distribution of time T + in x > when F < 0. Also 
T + + T~ = i hence a simple shift of the random variable 
yields T~ =t — T + , and hence also statistics for the case 
F < 0. 



D. Diffusion in an Unstable Force Field 

We consider a particle in an unstable force field 

3 :;s (^) 
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where F + > and F_ > 0. For this case the particle will 
eventually escape either to +oo or — oo, and the random 
walk is not recurrent. The survival probabilities in left 
and right domains are 

1 



where F + > and _F_ > 0. The random walk is recur- 
rent. Using an approach similar to one used in previous 
subsection 



W±(s) = - {1 - exp [ T F ± (1 + VI + 4s«r£) x / (2k b T)] } 

(75) 

where 



G s (u) 



s-\-u s 

i+ (s + u)+ (s) ' 



(83) 



(k b Tf 



Using Eq. we find 



G s (u) 



7+(s+") , 7-(s) 
s-\-u s 

7+ (s + u) + 7_ (s) 



(76) 



(77) 



where 



(78) 

Using the definition Eq. 0, the small s behavior of 
Eq. (|75l) gives the survival probabilities in the + and — 
domains, x > and x < respectively. In the limit of 
long times 



e±(x ) = 1 - exp ( =F 



_ F±x 
khT 



(79) 



where xq > for + and xq < for — . Hence according 
to the rather general Eqs. I|l(jll7fl we find 



G t (T+) 



F + 



F_ +F. 



-6 (T+ -t) + 



F++F- 



5{T+). (80) 



This long time behavior exhibits the same behavior for 
the normal diffusion a = 1 as for the anomalous case 
a < 1. Note that Eqs. 1)161171) where derived for normal 
diffusion however one can show that they are valid also 
for the anomalous diffusion case. To see this use the small 
s, u expansion of Eq. I|77|l which gives 



G s (u) 



1 



F+ + F- 



s 



(81) 



which is the double inverse Laplace transform of Eq. I|80|) 
and is independent of the parameters a, D a and T. 
For short times t <C min (r_, r+) we use the large s, u 

behavior of Eq. (JZZJ and find G t (T+) ~ \5 a/2 (l, 

Thus for short times the PDF of T + is the symmetric 
Lamperti PDF which is independent of all the parameters 
of the problem except for a. At the early stages of the 
dynamics the diffusion process not the drift is the most 
important, and hence forces are not relevant. 



E. Diffusion in Binding Force Field 

We consider a particle in a stable force field 
F(x)-- 



-F+ x>0 
F_ x<0 



(82) 



and £± (s) = F± (l - y/1 + 4s«r|) 
For small s.iiwe have 



G s (u) 



I . \a—l F_ ct — 1 

(S + U) jr; +S a 



and hence when t is large 



G t (t+) 



(84) 



(85) 



This is in agreement with our more general results 
Eqs. (HOI Ell). For short times we have G t {T + ) ~ 

\&a/2 flj which is similar to the behavior of the un- 
stable field discussed in previous sub-section. 



V. DISCUSSION 

Statistics of occupation times for binding external fields 
exhibits in the limit of long times an ergodic behavior 
when the diffusion is normal, or weak ergodicity break- 
ing Eqs. H40I44|) when diffusion is anomalous. We estab- 
lished a link between weak ergodicity breaking and frac- 
tional calculus. The exponent a in the fractional deriva- 
tive d a /dt a enters in Eq. 140|) describing the non-ergodic 
properties of the residence times. Since many processes 
and systems are modeled today using the fractional cal- 
culus approach, it is not out of the question that weak 
ergodicity breaking has many applications, and is wide 
spread. We can say that at-least one must treat with 
care, results obtained using fractional kinetic equations, 
since they describe only ensemble averages, not time av- 
erages. 

For binding external fields our results are in full agree- 
ment with those derived recently, by Bel and the au- 
thor |19|. There a continuous time random walk process 
was considered. Technically the methods used to treat 
the two problems arc different. For the fractional frame- 
work, a differential equation, Eq. (|2~4")l for the occupation 
times is derived and solved, which yields the weakly non 
ergodic properties of the system, while for the CTRW 
certain recursive relations must be solved 20]. Our work 
shows how the fractional framework, which is the contin- 
uum limit of the CTRW (and in this sense simpler) can 
be used to obtain statistics of residence times, and for 
binding fields weak ergodicity breaking. 

We found a general relation between the problem of 
occupation times and the problem of first passage times, 
Eq. 0. Mathematically the problem of first passage 
time is described in terms of a differential equation whose 



11 



solution depends on a single time t as a parameter, while 
the equation for residence time depends on two times 
T + and t, hence the approach based on calculation of 
first passage times instead of a direct calculation of dis- 
tribution of residence times seems to us useful. Besides, 
statistics of first passage times in external fields is a well 
investigated problem and this knowledge was used here. 
We showed that from Eq. Q general properties of res- 
idence times can be easily derived. Luckily the relation 
Eq. Q is valid for both the non-Markovian sub-diffusive 
and normal processes, and hence has some general va- 
lidity. From our derivation of this relation, we see that 
a key ingredient for this relation to be valid is the re- 
newal property of the underlying random walk. As far 
as we know statistics of residence times for more gen- 
eral non-Markovian and non renewal dynamics was not 
investigated so far [45f . 



As mentioned, for binding fields we confirm old re- 
sults in the limit of long measurement times. New re- 
sults for non-binding force fields and also for non asymp- 
totic times are also found. For example, for the simplest 
case of anomalous diffusion in free space, we found that 
the PDF of occupation times in half space is a symmet- 
ric Lamperti PDF with index a/2. This is the natural 
generalization of the well known arcsine law. Hence we 
showed that the fractional kinetic framework is indeed a 
natural generalization of the ordinary diffusion process. 



Other behaviors are found when diffusion is not free. 
For example for short times we expect rather generally 
that occupation times statistics is described by a sym- 
metric Lamperti PDF, since for short times drift by force 
fields is a slow process if compared with the diffusion 
process. For long times we get either (i) weak ergodicity 
breaking if the force field is binding, or (ii) occupation 
times described by Eq. (|16|) if both the left and right 
random walks are not recurrent (note that this equation 
is also valid for anomalous sub diffusion), or (iii) when 
random walk in the left and right domain is not bounded 
and recurrent, in such a way that the average first return 
time is infinite, then a Lamperti PDF describes statistics 
of occupation times in the limit of long times. For normal 
diffusion process the asymmetry parameter is determined 
by Eq. I|19|) . We also analyzed other cases where the left 
random walk is recurrent but not bounded and the right 
random walk is bounded, showing deviations from simple 
scaling. 
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VI. APPENDIX A 

The Fokker-Planck equation for a particle in a force 
field reads 



dc(x, t) 



Of 



DiL-ppc(x, t) 



with the operator 



Lpp — 



d 2 



d F(x) 



dx 2 dx k B T ' 



(86) 



(87) 



where D\ is the diffusion coefficient. To prepare for the 
FFPE we rewrite Eq. I|86[l in an integral form 



c(x, t) — S(x — Xq) = D\ oil Lppc(x, t) 



(88) 



where S(x — xq) are the initial conditions. According to 
the fractional kinetic approach, we must replace the inte- 
gral qI^ with a fractional Riemann-Liouville integration, 
defined as an operation on a function Z(t) according to 



I?Z(t) 



Z(t>) 



r(a)7 (t-t'Y 



-dt\ 



(89) 



and for our purpose < a < 1 . The FFPE in its integral 
form is 

c(x, t) - S(x - x ) = D a I?L FP c(x, t) (90) 

where D a is a generalized diffusion coefficient. Several 
authorspresent this equation in different ways. Some- 
times further differentiation with respect to time is 
made in Eq. I|9U|I to return to a fractional differential 
equation instead of the fractional integral form. How- 
ever, following work of Gorenflo and Mainardi e.g. |34| 
there is now growing use of Caputo symbols which are 
more elegant. Such symbols are used in Eq. (IL'-il) which 
has the same meaning as Eq. (|9Up. Later we use the 
Laplace t — + s transform of the FFPE equation (|90|l 



sc(x, s) — 5(x — xq) = DaS 1 a Lppc(x,s). 



(91) 



We see that the solution of FFPE with a < 1 is related to 
the solution of the Fokker-Planck equation when a — 1 
in Laplace space. To obtain the solution of the FFPE 
from the solution of the usual Fokker-Planck equation 
we must make a replacement D\ — > D a s 1 ~ a 17]. This 
similarity transformation in s space can be inverted to 
real time |35| . and with it one can obtain a solution of 
the FFPE once the corresponding solution of the usual 
Fokker-Planck equation is known. This transformation is 
related to subordination and the inverse Levy transform 

^: 

Similarly the survival probability for the normal dif- 
fusion case, and the fractional case are related to each 
other, by a simple transformation in Laplace space |35| . 
According to Eq. (57) in Ref. HU 



(92) 
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where W+ X0 (s) is the survival probability for the frac- 
tional a < 1 or the normal case a = 1. Using Eqs. JS] 



I92|l it is easy to prove the validity of Eq. 1251) . 
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